analysis_tag <- '01_jaeger_descriptive'

library(pheatmap)
library(knitr)
library(scater)
library(Seurat)
library(Cairo)  
library(biomaRt)
library(dplyr)
## library(reshape2)
#library(iSEE)
library(Matrix)
library(scran)
## library(mvoutlier)
## library(shiny)
## library(kableExtra)
## library(velocyto.R)
## library(cowplot)
## library(corrplot)
library(DT)
library(viridis)
## library("CellMixS")

ac <- function(col, alpha = 1) {
    
    apply(sapply(col, col2rgb)/255, 2, function(x) rgb(x[1], 
        x[2], x[3], alpha = alpha))
}

Data load

Raw counts retrieval. Assumes all matrices share gene id and gene names

fns <- list.files(file.path("..", "data"), "*gz")
ec <- list()

for (fn in fns) {
    tag <- gsub("_rawCount.txt.gz", "", basename(fn))
    ec[[tag]] <- read.table(file.path("..", "data", fn), sep = "\t", 
        header = TRUE)
    rownames(ec[[tag]]) <- sprintf("%s_%s", ec[[tag]]$gene_id, 
        ec[[tag]]$gene_name)
    gene_ids <- ec[[tag]]$gene_id
    ec[[tag]] <- ec[[tag]][, -c(1, 2)]
    
}

sapply(ec, function(x) head(colnames(x)))

     Ascl1_A_12wk     Ascl1_A_5d            Ascl1_B_12wk                 
[1,] "Ascl1_12wk_A01" "Ascl1_d4_plate1_A01" "Ascl1_12wk_2_20190514BJ_A01"
[2,] "Ascl1_12wk_A02" "Ascl1_d4_plate1_A02" "Ascl1_12wk_2_20190514BJ_A02"
[3,] "Ascl1_12wk_A03" "Ascl1_d4_plate1_A03" "Ascl1_12wk_2_20190514BJ_A03"
[4,] "Ascl1_12wk_A04" "Ascl1_d4_plate1_A04" "Ascl1_12wk_2_20190514BJ_A04"
[5,] "Ascl1_12wk_A05" "Ascl1_d4_plate1_A05" "Ascl1_12wk_2_20190514BJ_A05"
[6,] "Ascl1_12wk_A06" "Ascl1_d4_plate1_A06" "Ascl1_12wk_2_20190514BJ_A06"
     Ascl1_B_5d                   Gli1_A_12wk               
[1,] "Ascl1_5d_II_20190403BJ_A01" "Gli1_12wk_20180822BJ_A01"
[2,] "Ascl1_5d_II_20190403BJ_A02" "Gli1_12wk_20180822BJ_A02"
[3,] "Ascl1_5d_II_20190403BJ_A03" "Gli1_12wk_20180822BJ_A03"
[4,] "Ascl1_5d_II_20190403BJ_A04" "Gli1_12wk_20180822BJ_A04"
[5,] "Ascl1_5d_II_20190403BJ_A05" "Gli1_12wk_20180822BJ_A05"
[6,] "Ascl1_5d_II_20190403BJ_A06" "Gli1_12wk_20180822BJ_A06"
     Gli1_A_5d                     Gli1_B_5d              
[1,] "Gli1_20180503BJ_Plate_I_B13" "X20190219_Gli1_5d_A01"
[2,] "Gli1_20180503BJ_Plate_I_B14" "X20190219_Gli1_5d_A02"
[3,] "Gli1_20180503BJ_Plate_I_B15" "X20190219_Gli1_5d_A03"
[4,] "Gli1_20180503BJ_Plate_I_B16" "X20190219_Gli1_5d_A04"
[5,] "Gli1_20180503BJ_Plate_I_B18" "X20190219_Gli1_5d_A05"
[6,] "Gli1_20180503BJ_Plate_I_B20" "X20190219_Gli1_5d_A06"

Checking all rownames are equivalent

for (item in names(ec)) {
    for (item2 in names(ec)) {
        stopifnot(all(rownames(ec[[item]]) == rownames(ec[[item2]])))
    }
}

Metadata (currently just the plate)

meta <- list()
for (tag in names(ec)) {
    meta[[tag]] <- data.frame(cell = colnames(ec[[tag]]), plate = tag)
    rownames(meta[[tag]]) <- meta[[tag]]$cell
}

Annotation, preQC

sce <- list()

for (item in names(ec)) {
    stopifnot(all(rownames(meta[[item]]) == colnames(ec[[item]])))
    sce[[item]] <- SingleCellExperiment(list(counts = as.matrix(ec[[item]])), 
        colData = meta[[item]])
    
}

Checking all rownames are identical

for (item in names(sce)) {
    for (item2 in names(sce)) {
        stopifnot(all(rownames(sce[[item]]) == rownames(sce[[item2]])))
    }
}

Genes annotation

Getting gene feature annotation. Assumes all data matrices share the same ensembl GIDs, and that this is mouse data. Only queries it the first time.

if (!file.exists(file.path("..", "data", "gene_annotation.RData"))) {
    ## ids = gene_ids
    filters = "ensembl_gene_id"
    attributes = c(filters, "chromosome_name", "gene_biotype", 
        "start_position", "end_position")
    biomart = "ENSEMBL_MART_ENSEMBL"
    dataset = "mmusculus_gene_ensembl"
    host = "www.ensembl.org"
    
    bmart <- biomaRt::useMart(biomart = biomart, dataset = dataset, 
        host = host)
    
    feature_info <- biomaRt::getBM(attributes = attributes, filters = filters, 
        values = gene_ids, mart = bmart)
    
    mm <- match(gene_ids, feature_info[[filters]])
    feature_info_full <- feature_info[mm, ]
    save(feature_info_full, file = file.path("..", "data", "gene_annotation.RData"))
} else {
    load(file.path("..", "data", "gene_annotation.RData"))
    
}

for (item in names(sce)) {
    object <- sce[[item]]
    old_rdata <- rowData(object)
    keep_cols <- !(colnames(old_rdata) %in% colnames(feature_info_full))
    new_rdata <- cbind(old_rdata[, keep_cols], feature_info_full)
    rowData(object) <- new_rdata
    
    rowData(object)$ensembl_gene_id <- gene_ids
    sce[[item]] <- object
    rm(object, old_rdata, keep_cols, new_rdata)
}

Mitochondrial/spike ins

Looking for mitochondrial/spike ins for filtering out cells afterwards.

We exclude genes that are not expressed in any cell.

for (item in names(sce)) {
    
    sce[[item]] <- sce[[item]][rowSums(counts(sce[[item]]) > 
        0) > 0, ]
    
    is.mt <- grepl("MT", rowData(sce[[item]])$chromosome_name)
    is.spike <- grepl("ERCC-", rownames(rowData(sce[[item]])))
    isSpike(sce[[item]], "Spike") <- is.spike
    
    
    dim(sce[[item]])
    ## sce[[item]] <- normalize(sce[[item]])
    sce[[item]] <- calculateQCMetrics(sce[[item]], feature_controls = list(ERCC = is.spike, 
        mito = is.mt))
}

Overall quality metrics

Mind the bump in library sizes, are these cells duplets or different anyhow?

for (item in names(sce)) {
    print(item)
    par(mfrow = c(2, 2), mar = c(5.1, 4.1, 0.1, 0.1))
    hist(sce[[item]]$total_counts/1000, xlab = "Library sizes (thousands)", 
        main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
    hist(sce[[item]]$total_features_by_counts, xlab = "Number of expressed genes", 
        main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
    hist(sce[[item]]$pct_counts_Spike, xlab = "ERCC proportion (%)", 
        ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
    hist(sce[[item]]$pct_counts_mito, xlab = "Mitochondrial proportion (%)", 
        ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
}

[1] "Ascl1_A_12wk"

[1] "Ascl1_A_5d"

[1] "Ascl1_B_12wk"

[1] "Ascl1_B_5d"

[1] "Gli1_A_12wk"

[1] "Gli1_A_5d"

[1] "Gli1_B_5d"

for (item in names(sce)) {
    print(item)
    par(mfrow = c(2, 2), mar = c(5.1, 4.1, 1.1, 1.1))
    
    plot(density(sce[[item]]$total_counts/1000), xlab = "library sizes (thousands)", 
        main = "")
    
    rug(sce[[item]]$total_counts/1000)
    
    plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$pct_counts_mito, 
        pch = 20, ylab = "library size (thousands)", xlab = "mitochondrial proportion (%)", 
        col = ac("black", 0.5))
    
    plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$total_features_by_counts, 
        pch = 20, ylab = "library size (thousands)", xlab = "number of genes", 
        col = ac("black", 0.5))
    
    plot(y = sce[[item]]$pct_counts_mito, x = sce[[item]]$total_features_by_counts, 
        pch = 20, xlab = "number of genes", ylab = "mitochondrial proportion (%)", 
        col = ac("black", 0.5))
}

[1] "Ascl1_A_12wk"

[1] "Ascl1_A_5d"

[1] "Ascl1_B_12wk"

[1] "Ascl1_B_5d"

[1] "Gli1_A_12wk"

[1] "Gli1_A_5d"

[1] "Gli1_B_5d"

Outlier removal

Removal of outliers for the library size and the number of expressed features and mitochondrial proportions.

for (item in names(sce)) {
    print(item)
    
    libsize.drop <- isOutlier(sce[[item]]$total_counts, nmads = 3, 
        type = "lower", log = TRUE)
    ## feature.drop <-
    ## isOutlier(sce[[item]]$total_features_by_counts, nmads=3,
    ## type='lower', log=TRUE)
    feature.drop <- sce[[item]]$total_features_by_counts >= 5000 | 
        sce[[item]]$total_features_by_counts < 1500
    
    spike.drop <- isOutlier(sce[[item]]$pct_counts_Spike, nmads = 3, 
        type = "higher")
    mito.drop <- isOutlier(sce[[item]]$pct_counts_mito, nmads = 2, 
        type = "higher")
    
    manualmito.drop <- sce[[item]]$pct_counts_mito >= 5
    
    ## extra libsize outlier (manual) manuallibsize.drop <-
    ## sce[[item]]$total_counts < 100000
    
    sce[[item]] <- sce[[item]][, !(libsize.drop | feature.drop | 
        spike.drop | mito.drop | manualmito.drop)]
    print(data.frame(ByLibSize = sum(libsize.drop), ByFeature = sum(feature.drop), 
        BySpike = sum(spike.drop), ByMito = sum(mito.drop), ByManualMitoDrop = sum(manualmito.drop[!mito.drop]), 
        Remaining = ncol(sce[[item]])))
    
}

[1] "Ascl1_A_12wk"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1         1       309       0     19              269        60
[1] "Ascl1_A_5d"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1        24        31       0     43                0       308
[1] "Ascl1_B_12wk"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1        49        71       0     24                0       291
[1] "Ascl1_B_5d"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1        62        27       0     11                0       297
[1] "Gli1_A_12wk"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1        50       171       0     37                0       195
[1] "Gli1_A_5d"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1         8        26       0     41                8       204
[1] "Gli1_B_5d"
  ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1        44       215       0     67               54       155

Overall quality metrics after outlier removal

Cell QC (number of genes detected, number of reads, mitochondrial) after outlier removal.

for (item in names(sce)) {
    print(item)
    
    par(mfrow = c(2, 2), mar = c(5.1, 4.1, 0.1, 0.1))
    hist(sce[[item]]$total_counts/1000, xlab = "Library sizes (thousands)", 
        main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
    hist(sce[[item]]$total_features_by_counts, xlab = "Number of expressed genes", 
        main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
    hist(sce[[item]]$pct_counts_Spike, xlab = "ERCC proportion (%)", 
        ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
    hist(sce[[item]]$pct_counts_mito, xlab = "Mitochondrial proportion (%)", 
        ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
}

[1] "Ascl1_A_12wk"

[1] "Ascl1_A_5d"

[1] "Ascl1_B_12wk"

[1] "Ascl1_B_5d"

[1] "Gli1_A_12wk"

[1] "Gli1_A_5d"

[1] "Gli1_B_5d"

for (item in names(sce)) {
    print(item)
    
    par(mfrow = c(2, 2), mar = c(5.1, 4.1, 1.1, 1.1))
    
    plot(density(sce[[item]]$total_counts/1000), xlab = "library sizes (thousands)", 
        main = "")
    
    rug(sce[[item]]$total_counts/1000)
    
    plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$pct_counts_mito, 
        pch = 20, ylab = "library size (thousands)", xlab = "mitochondrial proportion (%)", 
        col = ac("black", 0.5))
    
    plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$total_features_by_counts, 
        pch = 20, ylab = "library size (thousands)", xlab = "number of genes", 
        col = ac("black", 0.5))
    
    plot(y = sce[[item]]$pct_counts_mito, x = sce[[item]]$total_features_by_counts, 
        pch = 20, xlab = "number of genes", ylab = "mitochondrial proportion (%)", 
        col = ac("black", 0.5))
}

[1] "Ascl1_A_12wk"

[1] "Ascl1_A_5d"

[1] "Ascl1_B_12wk"

[1] "Ascl1_B_5d"

[1] "Gli1_A_12wk"

[1] "Gli1_A_5d"

[1] "Gli1_B_5d"

QC plots

In this setting, feature controls are mitochondrial genes.

for (item in names(sce)) {
    print(item)
    
    print(plotHighestExprs(sce[[item]], exprs_values = "counts"))
    
    p1 <- plotColData(sce[[item]], x = "total_counts", y = "total_features_by_counts")
    p2 <- plotColData(sce[[item]], x = "pct_counts_feature_control", 
        y = "total_features_by_counts")
    p3 <- plotColData(sce[[item]], x = "pct_counts_feature_control", 
        y = "pct_counts_in_top_50_features")
    
    p4 <- plotColData(sce[[item]], x = "total_counts", y = "pct_counts_in_top_500_features")
    
    multiplot(p1, p2, cols = 2)
    multiplot(p3, p4, cols = 2)
}

[1] "Ascl1_A_12wk"

[1] "Ascl1_A_5d"

[1] "Ascl1_B_12wk"

[1] "Ascl1_B_5d"

[1] "Gli1_A_12wk"

[1] "Gli1_A_5d"

[1] "Gli1_B_5d"

Filter lowly expressed genes and normalize

for (item in names(sce)) {
    tryCatch({
        print(item)
        
        genes_keep <- rowSums(counts(sce[[item]]) > 1) >= 10
        sce[[item]] <- sce[[item]][genes_keep, ]
        sce[[item]] <- scran::computeSumFactors(sce[[item]])
        ## no spike ins here sce[[item]] <-
        ## scran::computeSpikeFactors(sce[[item]])
        sce[[item]] <- normalize(sce[[item]])
    }, error = function(x) print(x))
}

[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"

Cell cycle classification

Cell cycle phase estimation (mouse)

for (item in names(sce)) {
    print(item)
    
    mouse.pairs <- readRDS(system.file("exdata", "mouse_cycle_markers.rds", 
        package = "scran"))
    assignments <- cyclone(sce[[item]], mouse.pairs, gene.names = rowData(sce[[item]])$ensembl_gene_id)
    print(table(assignments$phase))
    
    colData(sce[[item]])$est_cell_cycle_phase <- assignments$phases
}

[1] "Ascl1_A_12wk"

 G1 G2M   S 
 65  50   8 
[1] "Ascl1_A_5d"

 G1 G2M   S 
276  27   6 
[1] "Ascl1_B_12wk"

 G1 G2M   S 
297   8   2 
[1] "Ascl1_B_5d"

 G1 G2M   S 
289  17   6 
[1] "Gli1_A_12wk"

 G1 G2M   S 
241  53   4 
[1] "Gli1_A_5d"

 G1 G2M   S 
201  11   7 
[1] "Gli1_B_5d"

 G1 G2M   S 
250  21  11

Dimensionality reduction visualizations

Unintegrated

Transformation to Seurat v3 objects. Please note that during data scaling we regress out both number of UMIs (in this case this means library size, because the protocol has no UMIs on it) and the mitochondrial proportions.

RANGE_PCS <- 1:30

## https://satijalab.org/seurat/pancreas_integration_label_transfer.html

so <- list()

for (item in names(sce)) {
    print(item)
    
    so[[item]] <- as.Seurat(sce[[item]])
    so[[item]]@meta.data$orig.ident <- so[[item]]@meta.data$plate
    Idents(so[[item]]) <- "plate"
}

[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"

First evaluating the unintegrated data structure

## so <- set.all.ident(so, id = identity)


merged <- merge(x = so[[1]], y = so[2:length(so)])
merged <- NormalizeData(object = merged)
merged <- FindVariableFeatures(object = merged)
merged <- ScaleData(object = merged)
merged <- RunPCA(object = merged)

PC_ 1 
Positive:  ENSMUSG00000006333-Rps9, ENSMUSG00000003657-Calb2, ENSMUSG00000042436-Mfap4, ENSMUSG00000037984-Neurod6, ENSMUSG00000023927-Satb1, ENSMUSG00000038642-Ctss, ENSMUSG00000024754-Tmem2, ENSMUSG00000036905-C1qb, ENSMUSG00000030707-Coro1a, ENSMUSG00000020135-Apc2 
       ENSMUSG00000026360-Rgs2, ENSMUSG00000025400-Tac2, ENSMUSG00000036887-C1qa, ENSMUSG00000036896-C1qc, ENSMUSG00000031822-Gse1, ENSMUSG00000028581-Laptm5, ENSMUSG00000058715-Fcer1g, ENSMUSG00000024621-Csf1r, ENSMUSG00000029771-Irf5, ENSMUSG00000030579-Tyrobp 
       ENSMUSG00000021876-Rnase4, ENSMUSG00000052336-Cx3cr1, ENSMUSG00000025511-Tspan4, ENSMUSG00000022619-Mapk8ip2, ENSMUSG00000031821-Gins2, ENSMUSG00000052160-Pld4, ENSMUSG00000074622-Mafb, ENSMUSG00000021423-Ly86, ENSMUSG00000031530-Dusp4, ENSMUSG00000007021-Syngr3 
Negative:  ENSMUSG00000026424-Gpr37l1, ENSMUSG00000020591-Ntsr2, ENSMUSG00000045092-S1pr1, ENSMUSG00000050953-Gja1, ENSMUSG00000030605-Mfge8, ENSMUSG00000006205-Htra1, ENSMUSG00000007097-Atp1a2, ENSMUSG00000022132-Cldn10, ENSMUSG00000023913-Pla2g7, ENSMUSG00000028517-Plpp3 
       ENSMUSG00000017390-Aldoc, ENSMUSG00000054252-Fgfr3, ENSMUSG00000029309-Sparcl1, ENSMUSG00000022037-Clu, ENSMUSG00000004892-Bcan, ENSMUSG00000029838-Ptn, ENSMUSG00000040055-Gjb6, ENSMUSG00000035805-Mlc1, ENSMUSG00000030428-Ttyh1, ENSMUSG00000041329-Atp1b2 
       ENSMUSG00000030495-Slc7a10, ENSMUSG00000055254-Ntrk2, ENSMUSG00000017009-Sdc4, ENSMUSG00000028128-F3, ENSMUSG00000024411-Aqp4, ENSMUSG00000032281-Acsbg1, ENSMUSG00000058135-Gstm1, ENSMUSG00000046240-Hepacam, ENSMUSG00000006522-Itih3, ENSMUSG00000005360-Slc1a3 
PC_ 2 
Positive:  ENSMUSG00000036896-C1qc, ENSMUSG00000036887-C1qa, ENSMUSG00000036905-C1qb, ENSMUSG00000028581-Laptm5, ENSMUSG00000052336-Cx3cr1, ENSMUSG00000030579-Tyrobp, ENSMUSG00000038642-Ctss, ENSMUSG00000054675-Tmem119, ENSMUSG00000021665-Hexb, ENSMUSG00000015852-Fcrls 
       ENSMUSG00000023992-Trem2, ENSMUSG00000048163-Selplg, ENSMUSG00000058715-Fcer1g, ENSMUSG00000059498-Fcgr3, ENSMUSG00000024621-Csf1r, ENSMUSG00000036908-Unc93b1, ENSMUSG00000036353-P2ry12, ENSMUSG00000020101-Vsir, ENSMUSG00000051504-Siglech, ENSMUSG00000040747-Cd53 
       ENSMUSG00000027848-Olfml3, ENSMUSG00000021423-Ly86, ENSMUSG00000052160-Pld4, ENSMUSG00000018774-Cd68, ENSMUSG00000002111-Spi1, ENSMUSG00000029919-Hpgds, ENSMUSG00000036362-P2ry13, ENSMUSG00000021666-Gfm2, ENSMUSG00000015396-Cd83, ENSMUSG00000040229-Gpr34 
Negative:  ENSMUSG00000003657-Calb2, ENSMUSG00000020135-Apc2, ENSMUSG00000068740-Celsr2, ENSMUSG00000048078-Tenm4, ENSMUSG00000042436-Mfap4, ENSMUSG00000047787-Flrt1, ENSMUSG00000032181-Scg3, ENSMUSG00000027309-4930402H24Rik, ENSMUSG00000025400-Tac2, ENSMUSG00000032086-Bace1 
       ENSMUSG00000031517-Gpm6a, ENSMUSG00000037984-Neurod6, ENSMUSG00000023927-Satb1, ENSMUSG00000025203-Scd2, ENSMUSG00000027282-Mtch2, ENSMUSG00000022619-Mapk8ip2, ENSMUSG00000046573-Lyrm4, ENSMUSG00000027692-Tnik, ENSMUSG00000058254-Tspan7, ENSMUSG00000019124-Scrn1 
       ENSMUSG00000024754-Tmem2, ENSMUSG00000031617-Tmem184c, ENSMUSG00000007021-Syngr3, ENSMUSG00000040111-Gramd1b, ENSMUSG00000006333-Rps9, ENSMUSG00000031821-Gins2, ENSMUSG00000031822-Gse1, ENSMUSG00000043733-Ptpn11, ENSMUSG00000033998-Kcnk1, ENSMUSG00000000202-Btbd17 
PC_ 3 
Positive:  ENSMUSG00000020423-Btg2, ENSMUSG00000052684-Jun, ENSMUSG00000006333-Rps9, ENSMUSG00000024190-Dusp1, ENSMUSG00000052837-Junb, ENSMUSG00000030737-Slco2b1, ENSMUSG00000031821-Gins2, ENSMUSG00000028976-Slc2a5, ENSMUSG00000036931-Nfkbid, ENSMUSG00000044786-Zfp36 
       ENSMUSG00000027399-Il1a, ENSMUSG00000031822-Gse1, ENSMUSG00000035356-Nfkbiz, ENSMUSG00000024073-Birc6, ENSMUSG00000026923-Notch1, ENSMUSG00000068740-Celsr2, ENSMUSG00000026628-Atf3, ENSMUSG00000000078-Klf6, ENSMUSG00000048078-Tenm4, ENSMUSG00000038393-Txnip 
       ENSMUSG00000019039-Dalrd3, ENSMUSG00000000028-Cdc45, ENSMUSG00000031613-Hpgd, ENSMUSG00000022817-Itgb5, ENSMUSG00000045730-Adrb2, ENSMUSG00000026965-Anapc2, ENSMUSG00000024401-Tnf, ENSMUSG00000002496-Tsc2, ENSMUSG00000027309-4930402H24Rik, ENSMUSG00000024975-Pdcd4 
Negative:  ENSMUSG00000092511-Gm20547, ENSMUSG00000090231-Cfb, ENSMUSG00000025498-Irf7, ENSMUSG00000074896-Ifit3, ENSMUSG00000017493-Igfbp4, ENSMUSG00000030560-Ctsc, ENSMUSG00000075014-Gm10800, ENSMUSG00000062488-Ifit3b, ENSMUSG00000035493-Tgfbi, ENSMUSG00000035385-Ccl2 
       ENSMUSG00000022037-Clu, ENSMUSG00000063903-Klk1, ENSMUSG00000033880-Lgals3bp, ENSMUSG00000060802-B2m, ENSMUSG00000038642-Ctss, ENSMUSG00000064373-Selenop, ENSMUSG00000033208-S100b, ENSMUSG00000002289-Angptl4, ENSMUSG00000090698-Apold1, ENSMUSG00000061232-H2-K1 
       ENSMUSG00000042659-Arrdc4, ENSMUSG00000036908-Unc93b1, ENSMUSG00000022257-Laptm4b, ENSMUSG00000071637-Cebpd, ENSMUSG00000036256-Igfbp7, ENSMUSG00000017390-Aldoc, ENSMUSG00000027199-Gatm, ENSMUSG00000036905-C1qb, ENSMUSG00000058715-Fcer1g, ENSMUSG00000023367-Tmem176a 
PC_ 4 
Positive:  ENSMUSG00000035686-Thrsp, ENSMUSG00000027004-Frzb, ENSMUSG00000027878-Notch2, ENSMUSG00000038668-Lpar1, ENSMUSG00000020932-Gfap, ENSMUSG00000037206-Islr, ENSMUSG00000052957-Gas1, ENSMUSG00000021250-Fos, ENSMUSG00000029413-Naaa, ENSMUSG00000015090-Ptgds 
       ENSMUSG00000059325-Hopx, ENSMUSG00000028195-Cyr61, ENSMUSG00000028214-Gem, ENSMUSG00000001025-S100a6, ENSMUSG00000067818-Myl9, ENSMUSG00000052684-Jun, ENSMUSG00000021702-Thbs4, ENSMUSG00000053398-Phgdh, ENSMUSG00000028583-Pdpn, ENSMUSG00000037031-Tspan15 
       ENSMUSG00000032348-Gsta4, ENSMUSG00000071637-Cebpd, ENSMUSG00000024190-Dusp1, ENSMUSG00000053519-Kcnip1, ENSMUSG00000028691-Prdx1, ENSMUSG00000059970-Hspa2, ENSMUSG00000023034-Nr4a1, ENSMUSG00000021186-Fbln5, ENSMUSG00000027253-Lrp4, ENSMUSG00000020423-Btg2 
Negative:  ENSMUSG00000033998-Kcnk1, ENSMUSG00000060961-Slc4a4, ENSMUSG00000038094-Atp13a4, ENSMUSG00000036949-Slc39a12, ENSMUSG00000022122-Ednrb, ENSMUSG00000055003-Lrtm2, ENSMUSG00000027200-Sema6d, ENSMUSG00000003657-Calb2, ENSMUSG00000040055-Gjb6, ENSMUSG00000020135-Apc2 
       ENSMUSG00000037348-Paqr7, ENSMUSG00000030495-Slc7a10, ENSMUSG00000026657-Frmd4a, ENSMUSG00000053004-Hrh1, ENSMUSG00000073424-Cyp4f15, ENSMUSG00000033849-B3galt2, ENSMUSG00000001260-Gabrg1, ENSMUSG00000037697-Ddhd1, ENSMUSG00000020734-Grin2c, ENSMUSG00000035376-Hacd2 
       ENSMUSG00000024411-Aqp4, ENSMUSG00000096054-Syne1, ENSMUSG00000002475-Abhd3, ENSMUSG00000022537-Tmem44, ENSMUSG00000010064-Slc38a3, ENSMUSG00000047787-Flrt1, ENSMUSG00000052387-Trpm3, ENSMUSG00000028766-Alpl, ENSMUSG00000000202-Btbd17, ENSMUSG00000032036-Kirrel3 
PC_ 5 
Positive:  ENSMUSG00000033720-Sfxn5, ENSMUSG00000058254-Tspan7, ENSMUSG00000017390-Aldoc, ENSMUSG00000005089-Slc1a2, ENSMUSG00000053398-Phgdh, ENSMUSG00000059325-Hopx, ENSMUSG00000020591-Ntsr2, ENSMUSG00000035805-Mlc1, ENSMUSG00000041556-Fbxo2, ENSMUSG00000032281-Acsbg1 
       ENSMUSG00000004902-Slc25a18, ENSMUSG00000028976-Slc2a5, ENSMUSG00000035686-Thrsp, ENSMUSG00000020932-Gfap, ENSMUSG00000089675-Ugt1a8, ENSMUSG00000090175-Ugt1a9, ENSMUSG00000089943-Ugt1a5, ENSMUSG00000090171-Ugt1a2, ENSMUSG00000090165-Ugt1a10, ENSMUSG00000089960-Ugt1a1 
       ENSMUSG00000090145-Ugt1a6b, ENSMUSG00000021792-Fam213a, ENSMUSG00000054545-Ugt1a6a, ENSMUSG00000026628-Atf3, ENSMUSG00000030737-Slco2b1, ENSMUSG00000026701-Prdx6, ENSMUSG00000027399-Il1a, ENSMUSG00000004558-Ndrg2, ENSMUSG00000045608-Dbx2, ENSMUSG00000030701-Plekhb1 
Negative:  ENSMUSG00000029718-Pcolce, ENSMUSG00000017344-Vtn, ENSMUSG00000025780-Itih5, ENSMUSG00000061353-Cxcl12, ENSMUSG00000042116-Vwa1, ENSMUSG00000029661-Col1a2, ENSMUSG00000036256-Igfbp7, ENSMUSG00000030717-Nupr1, ENSMUSG00000004665-Cnn2, ENSMUSG00000022548-Apod 
       ENSMUSG00000003617-Cp, ENSMUSG00000050212-Eva1b, ENSMUSG00000024620-Pdgfrb, ENSMUSG00000031538-Plat, ENSMUSG00000004951-Hspb1, ENSMUSG00000036814-Slc6a20a, ENSMUSG00000022090-Pdlim2, ENSMUSG00000031503-Col4a2, ENSMUSG00000027800-Tm4sf1, ENSMUSG00000063796-Slc22a8 
       ENSMUSG00000008999-Bmp7, ENSMUSG00000072941-Sod3, ENSMUSG00000017754-Pltp, ENSMUSG00000020467-Efemp1, ENSMUSG00000062980-Cped1, ENSMUSG00000053279-Aldh1a1, ENSMUSG00000024168-Tmem204, ENSMUSG00000001506-Col1a1, ENSMUSG00000089945-Pakap, ENSMUSG00000038729-Akap2

merged <- FindNeighbors(object = merged)
merged <- FindClusters(object = merged)

Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck

Number of nodes: 1862
Number of edges: 56945

Running Louvain algorithm...
Maximum modularity in 10 random starts: 0.8731
Number of communities: 14
Elapsed time: 0 seconds

merged <- RunTSNE(object = merged)
merged <- RunUMAP(object = merged, reduction = "pca", dims = RANGE_PCS)

p <- DimPlot(object = merged, reduction = "tsne", group.by = "plate")
print(p + theme(legend.position = "bottom"))

p2 <- DimPlot(object = merged, reduction = "umap", group.by = "plate", 
    label = TRUE, repel = TRUE) + NoLegend()

print(p2)

rm(merged)

Integrated

Second, CCA+MNN Seurat v3 integration

for (i in 1:length(so)) {
    so[[i]] <- NormalizeData(object = so[[i]], verbose = FALSE)
    so[[i]] <- FindVariableFeatures(object = so[[i]], selection.method = "vst", 
        nfeatures = 2000, verbose = FALSE)
}

jaeger.anchors <- FindIntegrationAnchors(object.list = so, dims = RANGE_PCS)

integrated <- IntegrateData(anchorset = jaeger.anchors, dims = RANGE_PCS)

DefaultAssay(object = integrated) <- "integrated"

## Run the standard workflow for visualization and clustering
integrated <- ScaleData(object = integrated, verbose = FALSE, 
    display.progress = FALSE, vars.to.regress = c("nUMI", "pct_counts_mito"))

integrated <- RunPCA(object = integrated, npcs = max(RANGE_PCS), 
    verbose = FALSE)
integrated <- RunTSNE(object = integrated, npcs = max(RANGE_PCS), 
    verbose = FALSE)

ElbowPlot(integrated)

RANGE_PCS_RED <- 1:10

integrated <- RunUMAP(object = integrated, reduction = "pca", 
    dims = RANGE_PCS_RED)

p1 <- DimPlot(object = integrated, reduction = "tsne", group.by = "plate")
p2 <- DimPlot(object = integrated, reduction = "umap", group.by = "plate", 
    label = TRUE, repel = TRUE) + NoLegend()
p3 <- DimPlot(object = integrated, reduction = "pca", group.by = "plate")

plot(p3)

plot(p1)

plot(p2)

Cell clustering

RES <- "integrated_snn_res.1"
res <- 1
integrated <- FindNeighbors(object = integrated, dims = RANGE_PCS_RED)

integrated <- FindClusters(object = integrated, resolution = res)

Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck

Number of nodes: 1862
Number of edges: 65447

Running Louvain algorithm...
Maximum modularity in 10 random starts: 0.7634
Number of communities: 11
Elapsed time: 0 seconds

integrated@active.assay <- "RNA"

Idents(object = integrated) <- RES

print(DimPlot(object = integrated, reduction = "umap", group.by = "plate"))

print(FeaturePlot(object = integrated, reduction = "umap", features = "pct_counts_mito"))

print(FeaturePlot(object = integrated, reduction = "umap", features = "nFeature_RNA"))

print(FeaturePlot(object = integrated, reduction = "umap", features = "nCount_RNA"))

print(DimPlot(object = integrated, reduction = "umap", group.by = RES))

markers <- FindAllMarkers(object = integrated, only.pos = TRUE, 
    min.pct = 0.25, thresh.use = 0.25)

Calculating cluster 0

Calculating cluster 1

Calculating cluster 2

Calculating cluster 3

Calculating cluster 4

Calculating cluster 5

Calculating cluster 6

Calculating cluster 7

Calculating cluster 8

Calculating cluster 9

Calculating cluster 10

DT::datatable(markers %>% as.data.frame() %>% dplyr::mutate_if(is.numeric, 
    funs(round(., 2))), extensions = c("Buttons", "FixedColumns"), 
    rownames = FALSE, options = list(dom = "Bfrtip", scrollX = TRUE, 
        fixedColumns = list(leftColumns = 1), buttons = c("csv", 
            "excel")))

markers_head <- as.data.frame(markers %>% group_by(cluster) %>% 
    top_n(n = 5, wt = avg_logFC))$gene

## graphics.off()
pheatmap(GetAssayData(object = integrated)[markers_head, ], color = viridis(100), 
    cluster_rows = FALSE, cluster_cols = TRUE, show_rownames = TRUE, 
    show_colnames = FALSE, clustering_distance_cols = "euclidean", 
    clustering_method = "ward.D2", clustering_distance_rows = "euclidean", 
    fontsize_row = 12, annotation_col = integrated@meta.data[, 
        c("plate", "nCount_RNA", "nFeature_RNA", "pct_counts_mito", 
            RES), drop = FALSE] %>% as.data.frame(), scale = "none")

Object export

save(sce, file = sprintf("%s_sce.RData", analysis_tag))
save(so, file = sprintf("%s_so.RData", analysis_tag))

Session

date()

[1] "Thu Jul 18 14:28:30 2019"

devtools::session_info()

─ Session info ──────────────────────────────────────────────────────────
 setting  value                       
 version  R version 3.6.0 (2019-04-26)
 os       Ubuntu 16.04.5 LTS          
 system   x86_64, linux-gnu           
 ui       X11                         
 language en_CA:en                    
 collate  en_CA.UTF-8                 
 ctype    en_CA.UTF-8                 
 tz       Europe/Zurich               
 date     2019-07-18                  

─ Packages ──────────────────────────────────────────────────────────────
 package              * version   date       lib source        
 AnnotationDbi          1.46.0    2019-05-02 [1] Bioconductor  
 ape                    5.3       2019-03-17 [1] CRAN (R 3.6.0)
 assertthat             0.2.1     2019-03-21 [1] CRAN (R 3.6.0)
 backports              1.1.4     2019-04-10 [1] CRAN (R 3.6.0)
 beeswarm               0.2.3     2016-04-25 [1] CRAN (R 3.6.0)
 bibtex                 0.4.2     2017-06-30 [1] CRAN (R 3.6.0)
 Biobase              * 2.44.0    2019-05-02 [1] Bioconductor  
 BiocGenerics         * 0.30.0    2019-05-02 [1] Bioconductor  
 BiocNeighbors          1.2.0     2019-05-02 [1] Bioconductor  
 BiocParallel         * 1.18.0    2019-05-03 [1] Bioconductor  
 BiocSingular           1.0.0     2019-05-02 [1] Bioconductor  
 biomaRt              * 2.40.0    2019-05-02 [1] Bioconductor  
 bit                    1.1-14    2018-05-29 [1] CRAN (R 3.6.0)
 bit64                  0.9-7     2017-05-08 [1] CRAN (R 3.6.0)
 bitops                 1.0-6     2013-08-17 [1] CRAN (R 3.6.0)
 blob                   1.1.1     2018-03-25 [1] CRAN (R 3.6.0)
 Cairo                * 1.5-10    2019-03-28 [1] CRAN (R 3.6.0)
 callr                  3.2.0     2019-03-15 [1] CRAN (R 3.6.0)
 caTools                1.17.1.2  2019-03-06 [1] CRAN (R 3.6.0)
 cli                    1.1.0     2019-03-19 [1] CRAN (R 3.6.0)
 cluster                2.0.9     2019-05-01 [1] CRAN (R 3.6.0)
 codetools              0.2-16    2018-12-24 [3] CRAN (R 3.6.0)
 colorspace             1.4-1     2019-03-18 [1] CRAN (R 3.6.0)
 colourpicker           1.0       2017-09-27 [1] CRAN (R 3.6.0)
 cowplot                0.9.4     2019-01-08 [1] CRAN (R 3.6.0)
 crayon                 1.3.4     2017-09-16 [1] CRAN (R 3.6.0)
 crosstalk              1.0.0     2016-12-21 [1] CRAN (R 3.6.0)
 data.table             1.12.2    2019-04-07 [1] CRAN (R 3.6.0)
 DBI                    1.0.0     2018-05-02 [1] CRAN (R 3.6.0)
 DelayedArray         * 0.10.0    2019-05-02 [1] Bioconductor  
 DelayedMatrixStats     1.6.0     2019-05-02 [1] Bioconductor  
 desc                   1.2.0     2018-05-01 [1] CRAN (R 3.6.0)
 devtools               2.0.2     2019-04-08 [1] CRAN (R 3.6.0)
 digest                 0.6.19    2019-05-20 [1] CRAN (R 3.6.0)
 dplyr                * 0.8.1     2019-05-14 [1] CRAN (R 3.6.0)
 dqrng                  0.2.1     2019-05-17 [1] CRAN (R 3.6.0)
 DT                   * 0.5       2018-11-05 [1] CRAN (R 3.6.0)
 dynamicTreeCut         1.63-1    2016-03-11 [1] CRAN (R 3.6.0)
 edgeR                  3.26.5    2019-06-21 [1] Bioconductor  
 evaluate               0.14      2019-05-28 [1] CRAN (R 3.6.0)
 fitdistrplus           1.0-14    2019-01-23 [1] CRAN (R 3.6.0)
 formatR                1.7       2019-06-11 [1] CRAN (R 3.6.0)
 fs                     1.3.0     2019-05-02 [1] CRAN (R 3.6.0)
 future                 1.12.0    2019-03-08 [1] CRAN (R 3.6.0)
 future.apply           1.2.0     2019-03-07 [1] CRAN (R 3.6.0)
 gbRd                   0.4-11    2012-10-01 [1] CRAN (R 3.6.0)
 gdata                  2.18.0    2017-06-06 [1] CRAN (R 3.6.0)
 GenomeInfoDb         * 1.20.0    2019-05-02 [1] Bioconductor  
 GenomeInfoDbData       1.2.1     2019-05-04 [1] Bioconductor  
 GenomicRanges        * 1.36.0    2019-05-02 [1] Bioconductor  
 ggbeeswarm             0.6.0     2017-08-07 [1] CRAN (R 3.6.0)
 ggplot2              * 3.2.0     2019-06-16 [1] CRAN (R 3.6.0)
 ggrepel                0.8.0     2018-05-09 [1] CRAN (R 3.6.0)
 ggridges               0.5.1     2018-09-27 [1] CRAN (R 3.6.0)
 globals                0.12.4    2018-10-11 [1] CRAN (R 3.6.0)
 glue                   1.3.1     2019-03-12 [1] CRAN (R 3.6.0)
 gplots                 3.0.1.1   2019-01-27 [1] CRAN (R 3.6.0)
 gridExtra              2.3       2017-09-09 [1] CRAN (R 3.6.0)
 gtable                 0.3.0     2019-03-25 [1] CRAN (R 3.6.0)
 gtools                 3.8.1     2018-06-26 [1] CRAN (R 3.6.0)
 hms                    0.4.2     2018-03-10 [1] CRAN (R 3.6.0)
 htmltools              0.3.6     2017-04-28 [1] CRAN (R 3.6.0)
 htmlwidgets            1.3       2018-09-30 [1] CRAN (R 3.6.0)
 httpuv                 1.5.1     2019-04-05 [1] CRAN (R 3.6.0)
 httr                   1.4.0     2018-12-11 [1] CRAN (R 3.6.0)
 ica                    1.0-2     2018-05-24 [1] CRAN (R 3.6.0)
 igraph                 1.2.4.1   2019-04-22 [1] CRAN (R 3.6.0)
 IRanges              * 2.18.1    2019-05-31 [1] Bioconductor  
 irlba                  2.3.3     2019-02-05 [1] CRAN (R 3.6.0)
 iSEE                 * 1.4.0     2019-05-02 [1] Bioconductor  
 jsonlite               1.6       2018-12-07 [1] CRAN (R 3.6.0)
 KernSmooth             2.23-15   2015-06-29 [3] CRAN (R 3.6.0)
 knitr                * 1.23      2019-05-18 [1] CRAN (R 3.6.0)
 labeling               0.3       2014-08-23 [1] CRAN (R 3.6.0)
 later                  0.8.0     2019-02-11 [1] CRAN (R 3.6.0)
 lattice                0.20-38   2018-11-04 [3] CRAN (R 3.6.0)
 lazyeval               0.2.2     2019-03-15 [1] CRAN (R 3.6.0)
 limma                  3.40.2    2019-05-17 [1] Bioconductor  
 listenv                0.7.0     2018-01-21 [1] CRAN (R 3.6.0)
 lmtest                 0.9-37    2019-04-30 [1] CRAN (R 3.6.0)
 locfit                 1.5-9.1   2013-04-20 [1] CRAN (R 3.6.0)
 lsei                   1.2-0     2017-10-23 [1] CRAN (R 3.6.0)
 magrittr               1.5       2014-11-22 [1] CRAN (R 3.6.0)
 MASS                   7.3-51.4  2019-03-31 [3] CRAN (R 3.6.0)
 Matrix               * 1.2-17    2019-03-22 [1] CRAN (R 3.6.0)
 matrixStats          * 0.54.0    2018-07-23 [1] CRAN (R 3.6.0)
 memoise                1.1.0     2017-04-21 [1] CRAN (R 3.6.0)
 metap                  1.1       2019-02-06 [1] CRAN (R 3.6.0)
 mgcv                   1.8-28    2019-03-21 [3] CRAN (R 3.6.0)
 mime                   0.7       2019-06-11 [1] CRAN (R 3.6.0)
 miniUI                 0.1.1.1   2018-05-18 [1] CRAN (R 3.6.0)
 munsell                0.5.0     2018-06-12 [1] CRAN (R 3.6.0)
 nlme                   3.1-139   2019-04-09 [3] CRAN (R 3.6.0)
 npsurv                 0.4-0     2017-10-14 [1] CRAN (R 3.6.0)
 pbapply                1.4-0     2019-02-05 [1] CRAN (R 3.6.0)
 pheatmap             * 1.0.12    2019-01-04 [1] CRAN (R 3.6.0)
 pillar                 1.4.1     2019-05-28 [1] CRAN (R 3.6.0)
 pkgbuild               1.0.3     2019-03-20 [1] CRAN (R 3.6.0)
 pkgconfig              2.0.2     2018-08-16 [1] CRAN (R 3.6.0)
 pkgload                1.0.2     2018-10-29 [1] CRAN (R 3.6.0)
 plotly                 4.9.0     2019-04-10 [1] CRAN (R 3.6.0)
 plyr                   1.8.4     2016-06-08 [1] CRAN (R 3.6.0)
 png                    0.1-7     2013-12-03 [1] CRAN (R 3.6.0)
 prettyunits            1.0.2     2015-07-13 [1] CRAN (R 3.6.0)
 processx               3.3.0     2019-03-10 [1] CRAN (R 3.6.0)
 progress               1.2.2     2019-05-16 [1] CRAN (R 3.6.0)
 promises               1.0.1     2018-04-13 [1] CRAN (R 3.6.0)
 ps                     1.3.0     2018-12-21 [1] CRAN (R 3.6.0)
 purrr                  0.3.2     2019-03-15 [1] CRAN (R 3.6.0)
 R.methodsS3            1.7.1     2016-02-16 [1] CRAN (R 3.6.0)
 R.oo                   1.22.0    2018-04-22 [1] CRAN (R 3.6.0)
 R.utils                2.8.0     2019-02-14 [1] CRAN (R 3.6.0)
 R6                     2.4.0     2019-02-14 [1] CRAN (R 3.6.0)
 RANN                   2.6.1     2019-01-08 [1] CRAN (R 3.6.0)
 RColorBrewer           1.1-2     2014-12-07 [1] CRAN (R 3.6.0)
 Rcpp                   1.0.1     2019-03-17 [1] CRAN (R 3.6.0)
 RCurl                  1.95-4.12 2019-03-04 [1] CRAN (R 3.6.0)
 Rdpack                 0.11-0    2019-04-14 [1] CRAN (R 3.6.0)
 remotes                2.0.4     2019-04-10 [1] CRAN (R 3.6.0)
 rentrez                1.2.2     2019-05-02 [1] CRAN (R 3.6.0)
 reshape2               1.4.3     2017-12-11 [1] CRAN (R 3.6.0)
 reticulate             1.12      2019-04-12 [1] CRAN (R 3.6.0)
 rintrojs               0.2.0     2017-07-04 [1] CRAN (R 3.6.0)
 rlang                  0.4.0     2019-06-25 [1] CRAN (R 3.6.0)
 rmarkdown              1.12      2019-03-14 [1] CRAN (R 3.6.0)
 ROCR                   1.0-7     2015-03-26 [1] CRAN (R 3.6.0)
 rprojroot              1.3-2     2018-01-03 [1] CRAN (R 3.6.0)
 RSQLite                2.1.1     2018-05-06 [1] CRAN (R 3.6.0)
 rsvd                   1.0.1     2019-06-02 [1] CRAN (R 3.6.0)
 Rtsne                  0.15      2018-11-10 [1] CRAN (R 3.6.0)
 S4Vectors            * 0.22.0    2019-05-02 [1] Bioconductor  
 scales                 1.0.0     2018-08-09 [1] CRAN (R 3.6.0)
 scater               * 1.12.2    2019-05-24 [1] Bioconductor  
 scran                * 1.12.1    2019-05-27 [1] Bioconductor  
 sctransform            0.2.0     2019-04-12 [1] CRAN (R 3.6.0)
 SDMTools               1.1-221.1 2019-04-18 [1] CRAN (R 3.6.0)
 sessioninfo            1.1.1     2018-11-05 [1] CRAN (R 3.6.0)
 Seurat               * 3.0.0     2019-04-15 [1] CRAN (R 3.6.0)
 shiny                  1.3.2     2019-04-22 [1] CRAN (R 3.6.0)
 shinyAce               0.3.3     2019-01-03 [1] CRAN (R 3.6.0)
 shinydashboard         0.7.1     2018-10-17 [1] CRAN (R 3.6.0)
 shinyjs                1.0       2018-01-08 [1] CRAN (R 3.6.0)
 SingleCellExperiment * 1.6.0     2019-05-02 [1] Bioconductor  
 statmod                1.4.32    2019-05-29 [1] CRAN (R 3.6.0)
 stringi                1.4.3     2019-03-12 [1] CRAN (R 3.6.0)
 stringr                1.4.0     2019-02-10 [1] CRAN (R 3.6.0)
 SummarizedExperiment * 1.14.0    2019-05-02 [1] Bioconductor  
 survival               2.44-1.1  2019-04-01 [3] CRAN (R 3.6.0)
 testthat               2.1.1     2019-04-23 [1] CRAN (R 3.6.0)
 tibble                 2.1.3     2019-06-06 [1] CRAN (R 3.6.0)
 tidyr                  0.8.3     2019-03-01 [1] CRAN (R 3.6.0)
 tidyselect             0.2.5     2018-10-11 [1] CRAN (R 3.6.0)
 tsne                   0.1-3     2016-07-15 [1] CRAN (R 3.6.0)
 usethis                1.5.0     2019-04-07 [1] CRAN (R 3.6.0)
 vipor                  0.4.5     2017-03-22 [1] CRAN (R 3.6.0)
 viridis              * 0.5.1     2018-03-29 [1] CRAN (R 3.6.0)
 viridisLite          * 0.3.0     2018-02-01 [1] CRAN (R 3.6.0)
 withr                  2.1.2     2018-03-15 [1] CRAN (R 3.6.0)
 xfun                   0.8       2019-06-25 [1] CRAN (R 3.6.0)
 XML                    3.98-1.20 2019-06-06 [1] CRAN (R 3.6.0)
 xtable                 1.8-4     2019-04-21 [1] CRAN (R 3.6.0)
 XVector                0.24.0    2019-05-02 [1] Bioconductor  
 yaml                   2.2.0     2018-07-25 [1] CRAN (R 3.6.0)
 zlibbioc               1.30.0    2019-05-02 [1] Bioconductor  
 zoo                    1.8-5     2019-03-21 [1] CRAN (R 3.6.0)

[1] /home/Shared/Rlib/release-3.9-lib
[2] /home/imallona/R/x86_64-pc-linux-gnu-library/3.6
[3] /usr/local/R/R-3.6.0/library

01 Baptiste Jaeger descriptive

Izaskun Mallona (Mark Robinson lab)